\section{Longitudinal Fixed Wing Flight Dynamics}


We are dealing here with the reduction of a fixed wing vehicle to the vertical plan.

\subsection{Frames}
The frames needed to solve the problem are depicted on figure \ref{fig: fdm_fw_longi_schematics}. We are going to use frames with $z$ pointing down in order to be consistent with notations that are commonly used for the 3D case.
\begin{figure}[h]
\centering
\fbox{
\includegraphics[width=\textwidth]{drawings/longitudinal_fixed_wing_schematics}
}
\caption{Frames used in the derivation of the dynamics.}
\label{fig: fdm_fw_longi_schematics}
\end{figure}
\begin{itemize}
  \item $R$ is the \emph{earth-fixed} frame, considered inertial. $x$ axis is pointing right, $z$ axis is pointing down.
  \item $B$ is the \emph{body-fixed} frame. It is located at the center of gravity of the vehicle and its $x$ axis is pointing toward the nose of the vehicle ($z$ is pointing down).
  \item $S$ is the \emph{stability} (or aerodynamic) frame. it is located at the center of mass of the vehicle and its $x$ axis is pointing in the direction of the airspeed of the vehicle ($z$ pointing down with respect to that).
  \item $P$ is the \emph{propulsion} frame. It is also a vehicle-fixed frame, but the $x$ axis is pointing in the direction of the axis of the motors.
\end{itemize}

Wind and its time derivative are noted 
\begin{equation}
  \vect{W}^R = \begin{pmatrix}W_x^R\\W_z^R\end{pmatrix} \quad \dot{\vect{W}}^R = \begin{pmatrix}\dot{W}_x^R\\\dot{W}_z^R\end{pmatrix}
\end{equation}


The transformation matrices between those frames are as follow:
\begin{equation}
  R^{R \rightarrow B} = \begin{pmatrix}\cos{\alpha+\gamma_a}&-\sin{\alpha+\gamma_a}\\\sin{\alpha+\gamma_a}&\cos{\alpha+\gamma_a}\end{pmatrix} \quad
  R^{R \rightarrow S} = \begin{pmatrix}\cos{\gamma_a}&-\sin{\gamma_a}\\\sin{\gamma_a}&\cos{\gamma_a}\end{pmatrix} \quad
  R^{S \rightarrow B} = \begin{pmatrix}\cos{\alpha}&-\sin{\alpha}\\\sin{\alpha}&\cos{\alpha}\end{pmatrix}
\end{equation}

\subsection{Forces}
%
\subsubsection*{Weight}
Weight is applying at the center of gravity of the vehicle and pointing toward $z_R$
\begin{equation}
  \vect{F}_{\text{weight}}^R = \begin{pmatrix}0\\m.g\end{pmatrix}
\end{equation}
%
\subsubsection*{Aerodynamic forces}
Aerodynamic forces are applying at the center of gravity of the vehicle and expressed in $S$ frame as follow:
\begin{equation}
  \vect{F}_{\text{aero}}^S = \begin{pmatrix}-D\\-L\end{pmatrix}
\end{equation}
Noting
\begin{equation}
  C_L = C_{L0} + C_{L\alpha}.\alpha + C_{Lq}.q + C_{Lde}.d_e
\end{equation}
lift and drag are expressed as
\begin{align}
  L &= \frac{1}{2}.\rho.v_a^2.S_\text{ref}.C_L\\
  D &= \frac{1}{2}.\rho.v_a^2.S_\text{ref}.\left(C_{D0} + C_{Dk1}.C_L + C_{Dk2}.C_L^2\right)
\end{align}
%
\subsubsection*{Propulsion forces}
Propulsion forces are applying at the location of the engine and expressed in $P$ frame as follow:
\begin{equation}
  \vect{F}_{\text{prop}}^P = \begin{pmatrix}T\\0\end{pmatrix}
\end{equation}
where
\begin{equation}
  T = T_\text{max}.(\frac{\rho}{\rho_i})^{n\rho}.(\frac{v_a}{v_{ai}})^{nv}.d_{th}
\end{equation}
%
\subsubsection*{Total}
The sum of all forces expressed in $B$ frame is computed as follow:
\begin{equation}
  \vect{F}^B = 
  R^{R \rightarrow B}.\begin{pmatrix}0\\m.g\end{pmatrix} +
  R^{S \rightarrow B}.\begin{pmatrix}-D\\-L\end{pmatrix} +
  R^{P \rightarrow B}.\begin{pmatrix}T\\0\end{pmatrix}
\end{equation}
%
\subsection{Momentum}
\begin{equation}
M^B = \frac{1}{2}.\rho.v_a^2.S_\text{ref}.c_\text{ref} \left(C_{m0} + C_{m\alpha}.\alpha + C_{mde}.d_e\right)
\end{equation}
\subsection{State Space Representation}
State Vector
\begin{equation}
  \vect{X} = \transp{\begin{pmatrix}x&z&v_a&\gamma_a&\alpha&q\end{pmatrix}}
\end{equation}

\begin{itemize}
  \item $x$ and $z$ are the coordinates of the center of gravity of the vehicle (in m)
  \item $v_a$ is the (norm of the) airspeed (in m/s)
  \item $\gamma_a$ is the air flight path (in rad)
  \item $\alpha$ is the angle of attack (in rad)
  \item $q$ is the angular velocity (in rad/s)
\end{itemize}


Input Vector
\begin{equation}
  \vect{U} = \transp{\begin{pmatrix}d_e&d_{th}\end{pmatrix}}
\end{equation}

\begin{itemize}
  \item $d_e$ is the deflection of the elevator (in rad)
  \item $d_{th}$ is position of the throttle lever (between 0 and 1)
\end{itemize}


\subsubsection*{Translational Dynamics}
The first two lines of the state space representation are obtained by applying kinematic relations;
\begin{align}
  \begin{pmatrix}\dot{x}\\\dot{z}\end{pmatrix} & = R^{S \rightarrow R}.\begin{pmatrix}va\\0\end{pmatrix} + \vect{W}^R\\
  \label{eq: longi_fw_state_space_xz}
  \begin{pmatrix}\dot{x}\\\dot{z}\end{pmatrix} & =
  \begin{pmatrix}
    v_a.\cos{\gamma_a} + W_x^R\\ -v_a.\sin{\gamma_a} + W_z^R
  \end{pmatrix}
\end{align}
%
The next two lines are obtained using Newton's second law for forces in earth frame
\begin{equation}\label{eq: longi_fw_newton_second_law_earth}
  \vect{F}^R = m.\dot{\vect{v}}^R
\end{equation}
%
The norm of the airspeed is expressed in terms of the components of the airspeed expressed in earth frame as follow:
\begin{equation}
  v_a = \sqrt{(v_{ax}^R)^2+(v_{az}^R)^2}
\end{equation}
Its time derivative is computed as follow:
\begin{equation}\label{eq: longi_fw_state_space_va}
  \dot{v}_a = \frac{v_{ax}^R.\dot{v}_{ax}^R + v_{az}^R.\dot{v}_{az}^R}{\sqrt{(v_{ax}^R)^2+(v_{az}^R)^2}} =
  \frac{v_{ax}^R.\dot{v}_{ax}^R + v_{az}^R.\dot{v}_{az}^R}{v_a}
\end{equation}
%
The air flight path angle $\gamma_a$ is expressed in terms of the components of the airspeed expressed in earth frame as follow:
\begin{equation}
  \gamma_a = -atan(\frac{v_{az}^R}{v_{ax}^R})
\end{equation}
Its time derivative is computed as:
\begin{equation}\label{eq: longi_fw_state_space_gamma}
  \dot{\gamma}_a = \frac{-v_{ax}^R.\dot{v}_{az}^R + v_{az}^R.\dot{v}_{ax}^R }{(v_{ax}^R)^2+(v_{az}^R)^2} = 
  \frac{-v_{ax}^R.\dot{v}_{az}^R + v_{az}^R.\dot{v}_{ax}^R }{v_{a}^2}
\end{equation}
%
The componnents of airspeed in earth frame are expressed in terms of inertial speed as
\begin{equation}
  \vect{v}_a^R = \vect{v}^R - \vect{W}^R
\end{equation}
The time derivative is then computed as:
\begin{equation}
  \dot{\vect{v}_a^R} = \dot{\vect{v}}^R - \dot{\vect{W}}^R
\end{equation}
Introducing equation \eqref{eq: longi_fw_newton_second_law_earth} leads to
\begin{equation}\label{eq: longi_fw_state_space_vad}
  \dot{\vect{v}_a^R} = \frac{1}{m}\vect{F}^R - \dot{\vect{W}}^R
\end{equation}


\subsubsection*{Rotational Dynamics}

The fifth line of the state space representation is obtained using kinematics:
\begin{equation}\label{eq: longi_fw_state_space_alpha}
  \dot{\alpha} = q - \dot{\gamma_a}
\end{equation}
The last line of the state space representation is obtained by applying Newton's second law for moments in body frame
\begin{equation}\label{eq: longi_fw_state_space_q}
  I.\dot{q} = \vect{M}^B
\end{equation}

\subsubsection*{Complete Dynamics}
Regrouping equations \ref{eq: longi_fw_state_space_xz}, \ref{eq: longi_fw_state_space_va} , \ref{eq: longi_fw_state_space_gamma}, \ref{eq: longi_fw_state_space_vad}, \ref{eq: longi_fw_state_space_alpha} and \ref{eq: longi_fw_state_space_q} leads to

\begin{align}
  \dot{\vect{X}} & = f(\vect{X}, \vect{U})\\
  \begin{pmatrix}\dot{x}\\\dot{z}\\\dot{v}_a\\\dot{\gamma}_a\\\dot{\alpha}\\\dot{q}\end{pmatrix} & = 
  \begin{pmatrix}
    v_a.\cos{\gamma_a}+W_x^R\\
    -v_a.\sin{\gamma_a}+W_z^R\\
    \frac{T}{m}.\cos{\alpha} -\frac{D}{m} -g.\sin{\gamma_a} - \dot{W}_x^R.\cos{\gamma_a} + \dot{W}_z^R.\sin{\gamma_a}\\
    \frac{1}{v_a}\left( \frac{T}{m}.\sin{\alpha} +\frac{L}{m} -g.\cos{\gamma_a} + \dot{W}_x^R.\sin{\gamma_a} + \dot{W}_z^R.\cos{\gamma_a}\right)\\
    q - \dot{\gamma}_a\\
    \frac{1}{I}.M^B
  \end{pmatrix}
\end{align}
